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1. Introduction 

The study of many-body properties of ultra-cold and 
low density atomic samples is a topic of considerable 
experimental and theoretical interest. Substantial 
advances in this new subfield of atomic, molecular, and 
optical physics have been reported in the recent past, the 
most spectacular one being without doubt the observa- 
tion of Bose-Einstein condensation in weakly interacting 
alkali vapors [1-3]. On the theoretical front, quantum 
field theories of the interaction between atoms and light 
fields are being developed [4-6], and are finding appli- 
cations in the analysis of detection schemes for Bose 
condensates, nonlinear atom optics, etc. 

While in the experiments reported so far, Bose- 
Einstein condensation has been detected via the velocity 
distribution of the atoms in the condensate, it is of inter- 
est to examine alternate techniques, such as optical 
methods, which could give clear and unambiguous 
signatures of the distinctions between a condensate and 
an incoherent atomic sample of the same density [7-12]. 
In this paper, we use the analogy between nonlinear 



optics and nonlinear atom optics to propose and analyze 
another measurement scheme sensitive to the coherence 
properties of the atomic sample under consideration. 
Based on this analogy, we have predicted in earlier work 
[13] that it should be possible to achieve atomic phase 
conjugation from an atomic condensate. Since phase 
conjugation is a coherent optical effect, which relies on 
the mutual coherence between the optical fields in- 
volved, we argue and demonstrate numerically that 
atomic phase conjugation is rapidly destroyed by 
incoherent processes such as spontaneous emission. 
Hence, it is not expected to occur for an incoherent 
atomic sample, and is a sensitive probe of the coherence 
properties and of the lifetime of the condensate. 

Section 2 discusses our physical model, and proposes 
an experimental scheme where a Bose condensate is 
diffracted from an optical grating for atoms to produce 
two counter propagating condensate components, as 
well as a weak ' 'probe" which is scattered off the 
condensate. In a manner reminiscent of the optical case, 
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this configuration can lead to matter-wave phase conju- 
gation via four-wave mixing. The basic approximations 
underlying our analysis, as well as the nonlinear master 
equation which models our system, are presented. 

Section 3 reviews the coherent dynamics of the atoms 
and demonstrates that indeed, phase conjugation is 
possible in principle. Because the atoms under consider- 
ation are modeled as two-level systems, though, the 
theory of phase conjugation is somewhat more compli- 
cated than in the situation of degenerate four-wave 
mixing familiar in nonlinear optics, and is more remi- 
niscent of the case where the polarization of the optical 
fields must be taken into account. Spontaneous emission 
is included in Sec. 4. In this case, the nonlinear 
Schrodinger equation used to describe the coherent 
effects is replaced by a nonlinear Hartree-Fock master 
equation, which is solved numerically. Spontaneous 
emission is shown to drastically reduce the amount of 
phase conjugation available, thereby demonstrating its 
sensitivity to the coherence of the condensate. In addi- 
tion, we compare these results to the solution of a 
Hartree nonlinear Monte Carlo wave functions numeri- 
cal scheme that we have recently proposed, thereby 
gaining useful intuition about the way a condensate's 
description departs from Hartree as dissipation is 
increased. Finally, Sec. 5 is a summary and conclusion. 
Calculational details are relegated to appendices. 

In order to avoid any confusion, we note now that the 
nonlinear Schrodinger equation that we use in this paper 
is not the Gross -Pi taevskii nonlinear Schrodinger 
equation familiar in the description of Bose conden- 
sates, but rather the nonlinear Schrodinger equation of 
atom optics [4,14]. In the former equation, the effective 
nonlinearity results from short-range interactions 
between ground state atoms, while in the latter case it 
results from the long-range, near-resonant dipole-dipole 
interaction between ground and excited atoms. Because 
the diffraction grating involved in our scheme induces 
real transitions between ground and excited atoms, this 
latter interaction is dominant, and hence we neglect the 
short range ground-ground interactions here. After 
having neglected this interaction, we do treat the dipole- 
dipole interaction as local, so that our nonlinear 
Schrodinger equation is itself local, just like the 
Gross-Pitaevskii equation. However, the physical 
situation it describes is quite different, and indeed, its 
Schrodinger field is vectorial rather than scalar. 



2. Physical Model 

In analogy to the optical case of phase conjugation by 
four- wave mixing, we aim to realize a situation where a 
weak probe wave interacts nonlinearly with a strong 



* 'pump' ' wave to establish a grating from which a coun- 
ter propagating pump wave is scattered to produce a 
phase-conjugate wave. A possible way to establish this 
geometry is shown in Fig. 1, which illustrates how a 
Bose condensate is diffracted by a matter waves grating 
produced by a periodic evanescent electromagnetic 
wave. 




Fig. 1. Configuration of the atomic beams incident on and diffracted 
from the evanescent light grating. 

If the condensate impinges on the grating at the 
Bragg (Littrow) angle, the incident beam will be back- 
diffracted into two components, one propagating in the 
direction opposite to the incident wave (first-order 
diffraction) and the other, resulting from zeroth-order 
diffraction, at some angle from these ''beams." We 
already note at this point several important differences 
between this situation and the optical case: First, the 
waves are now de Broglie matter waves, instead of elec- 
tromagnetic fields; second, these waves have several 
components, corresponding to all electronic states of 
relevance for the system at hand, so that the fields under 
consideration are vector fields; finally, due to two-body 
interactions between the atoms, their evolution is in- 
trinsically nonlinear, so that no additional nonlinear 
"crystal" is required. 

In general, atoms impinging on the grating in their 
ground electronic state will be diffracted in a superposi- 
tion of ground and excited states, with two major conse- 
quences: the spontaneous decay of the excited atoms is 
an incoherent process that destroys the coherence 
between the various interacting waves, a mechanism that 
is expected to reduce the amount of phase conjugation 
achievable. In addition, a spontaneously emitted photon 
can be reabsorbed by another atom in the sample. This 
is the physical mechanism leading to the near-resonant 
dipole-dipole interaction, which is the major source of 
nonlinearity in the gedanken experiment at hand. 
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Note that in this equation, we have dropped the imagi- 
nary potential V c . While both V and V c find their 
physical origin in spontaneous emission, V conserves the 
purity of the state, while in general V c does not conserve 
Tr(p 2 ) [14]. As such, it should really be considered a 
dissipative term, and consistently, be neglected in the 
coherent regime. 

In complete analogy with conventional nonlinear 
optics [21-23] we proceed by introducing slowly vary- 
ing forward F(r,t) = (F e (r, t), F g (r, 1)) and backward 
B(r,t) = (B e (r,tX B g (r,t)) propagating components 
of the effective single particle wave function <fi(r,t) 



<$>{r, t) = F(r, Oe** 3 "" + B (r, t)e 



-i(kz+o)t) 



(7) 



with \V 2 F e (r 9 t)\ «k\V F € (r,t)\ «k 2 \F € (r,t)\ and 
\V 2 B € (r 9 t)\ «k\VB € (r,t)\ «k 2 \B € (r,t)\, where 
£= [e, g } labels the electronic state of the atoms. Here, 
we assume that the main direction of propagation of the 
atoms is along the z-axis, see Fig. 1. When substituted 
into the nonlinear Schrodinger equation, Eq. (6) with the 
explicit form of the dipole-dipole potential V, this 
decomposition yields a system of four coupled paraxial 
wave equations 



jL + 2 k £)F e = iVlF e 



-iV (\F g \ 2 F e + \B g \ 2 F e + F g B*B e ) 



± + 2k-?- 



-iV (\F g \ 2 B e +\B g \ 2 B e + B g F*F e ) 



j. t + 2kg)F, = iVlF. 



-iV (\F e \ 2 F g +\B e \ 2 F g + F e B*B g ) 



j- t + 2k£)B g = iVlB g 



-iV Q {\F e \ 2 B g +\B e \ 2 B g + B e F*F g ), 



(8) 



where we have introduced the dimensionless wave- 
number k —> k/k , time t —> (o rec t, wavef unction 
(f>(y, z, t) -> cf)(y, z, t)/k and position r -> k r 
variables, as well as a dimensionless potential 
V — » V /(o rec . Equations (8) resemble the usual model 
equations for beam propagation in nonlinear optics [21], 



except that due to the vectorial character of the slowly 
varying amplitudes F(r,t) and B(r,t) and the 
exchange of excitation involved in the dipole-dipole 
interaction, the nonlinear term only couples partial 
waves corresponding to different electronic states. As 
such, these equations are reminiscent of those involving 
polarization coupling in nonlinear optics, although in 
the latter case, terms involving self-phase modulation 
would also appear. Still, the analogy with nonlinear 
optics is sufficient in that we can readily adapt the beam 
propagation method common in optics to the problem at 
hand. 

For numerical purposes, we restrict our attention to 
two dimensions and to the portion of the condensate 
confined by the planes z = and z = L and assume that 
the boundary conditions are provided by the waves 
incident on these planes from the outside of that 
volume. 1 If we assume that the Bose condensate consists 
of ground state atoms sufficiently far from the grating, 
we have 

F e (y, 0; t) - ^ e0 (0, t) + #„(0, t)e"* 

F g (y 9 0; t) - ^ (0, t) + & gs (0 9 t)c iKy 

B e (y,L;t) = 

B g (y 9 L;t) = & g0 (L,t) . (9) 



Here, we have introduced a spatial Fourier decomposi- 
tion of the transverse variation of the amplitudes to 
account for diffraction by the grating, which is assumed 
to leave the atoms in a superposition of their excited and 
ground electronic states. Hence £F e0 (z, t), ^ g o(z, t) and 
£ft g o(z, t) are the amplitudes of the forward and back- 
ward propagating pump beams along the z-axis, while 
f €S (z, t) and ^ s (z, t) are the amplitudes of the forward 
propagating probe waves propagating at a small angle to 
the z-axis such that K«k. 

For pump waves strong enough that the undepleted 
pump approximation holds, i.e., %>(z, t) — const, 
S8fl)(z, t) — const, this problem can be solved analyti- 
cally [21,22]. It was shown in Ref. [13] that the intrinsic 
nonlinearity of the system results in the ground state 
probe being associated to an excited state conjugate 
wave propagating with transverse wave number -k , and 
similarly, the excited probe leads to a ground state 
conjugate. 



The x -dimension is eliminated by expressing the wave function 
(f>t(r,t) as <f>(x)(f>{(y,z,t) where <f>(x) is taken to be a solution 
of a linear Schrodinger equation. As a result V Q is replaced by 
V J dx\<Kx)\ 4 Iko. 
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We have mentioned in the Introduction that the 
dynamics of a Bose condensate is described by the 
Gross-Pitaevskii nonlinear Schrodinger equation [15], 
in which the nonlinearity results from short-range 
two-body collisions with a scattering length of the order 
of the Bohr radius. In the presence of light fields, 
however, the two-body potential is dominated at low 
densities by the dipole-dipole interaction between 
ground and excited atoms. Neglecting the short-range 
interactions, subsequently approximating the dipole- 
dipole potential by a ' 'contact" potential, and perform- 
ing the Hartree approximation yields an effective single- 
particle nonlinear Schrodinger equation which forms 
the basis of nonlinear atom optics [4]. Although this 
equation has the same structure as the Gross-Pitaevskii 
equation, it describes substantially different physics. 

However, since spontaneous emission cannot be 
ignored when excited atomic levels are involved (and in 
addition is the cause of the dipole-dipole interaction) the 
nonlinear Schrodinger equation does not give an appro- 
priate description of the system. Rather, it must be re- 
placed by a master equation approach. This last equation 
results from the common procedure of eliminating the 
degrees of freedom of the continuum of modes of the 
electromagnetic field in the Born-Markoff approxima- 
tion. Invoking in addition the Hartree-Fock factorization 
Ansatz results in the effective single particle nonlinear 
master equation [16,14] 



H = — + ficoocr+cr-, 

1m 



(2) 



where p is the atomic center of mass momentum, m is 
the single particle mass, a) = k c is the atomic transi- 
tion frequency and a+ are the usual pseudospin raising 
and lowering operators. The dipole-dipole potential V 
takes the explicit form [17-19] 



V=4f (AM){yo(*o|r, 2 |)-i {l-3(c(r 12 ) • (e(d)) 2 } 



X y 2 (k \r l2 \)}((r + ® o-+ o-® <r + ) 
and the imaginary potential V c is 



(3) 



V c = i^f(N-l){jo(k \r l2 \)-^{l-3(e(r l2 y(e(d)) 2 } 



XJ2(^o|ri 2 |)}{o-®o- + - o- + ®o-}. 



(4) 



In this expression, j t , y t are modified spherical Bessel 
functions of the first and second kind respectively [20], 
N notes the number of atoms in the sample, 
r= \d 1 2 kl 1 3h is the spontaneous emission rate and d 
the dipole matrix element of the atomic transition. The 
relative position of the atoms is denoted r n and e (x) is 
a unit vector in the direction of the vector x. In the 
following, we approximate these potentials by a contact 
potential of the form 



^^ = -Ua\[H,p]\b) + (a\L nc p\b) 



dt 



- l - h J"d{i}[<a:l|(V+V c )|2:3>(p(2;i)p(3;l) 

+ p(3;b)p(2;l)) 

-(p(a;l)p(3;2)) + p(a;2)p(3;l)) 

<1:2|(V+V c )|i:3>]. (1) 



V+ V c = -^ 8 (r 12 )(<7 + ® a- + o- ® <r + ) 



iV c0 
+ —78 (r 12 )(<7- ® <7 + + <7 + ® <7") , (5) 



an approximation adequate for ultra-cold atoms when 
the thermal de Broglie wavelength of the atoms is much 
larger than an optical wavelength. 



Here p(a;b) are the matrix elements of the single 
particle density operator p, the "integration over 
numbers" such as d{i} means a summation over a 
complete set of quantum numbers for the single particle 
system, H is the single atom Hamiltonian, L nc p is the 
usual single atom spontaneous emission term, the 
two-body potential V describes the dipole-dipole inter- 
action, and V c is an imaginary potential which results 
from the elimination of the modes of the thermal elec- 
tromagnetic field in the Born-Markoff approximation. 
We specialize to the case of a system of two-level 
atoms with lower electronic level |g) and upper elec- 
tronic level \e). The single atom Hamiltonian then takes 
the form 



3. Coherent Regime 

Having established our model, we briefly review mat- 
ter waves phase conjugation in the coherent regime. In 
this case, the dynamics of the system is described by the 
Hartree nonlinear Schrodinger equation [4] 

m d&£) = ( d2 ( a \H\2)cf>(2) 

+ J dld2d3<a,l|V|2,3)0*(l)0(2)0(3) . (6) 
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We performed numerical simulations of the coupled 
paraxial wave equations, Eq. (8) for various boundary 
conditions involving different fractions r\ of excited 
atoms in the diffracted wave at z = 0. As a numerical 
check, we verified that the atomic flux through the z = 
and z = L planes is properly conserved. Fig. 2 shows an 
example of the build-up of the ground state phase conju- 
gate wave of amplitude £$ g (- k, 0, t), and the associated 
depletion of the pump amplitude ^(0, 0, t). 
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Fig. 2. Time- dependence of the backward propagating atomic 
beams| %(- k, 0, t)\ 2 — dashed line, and | %(0, 0, t)\ 2 — solid line, 
for k = k , k= lOito, L = 40/* , Vo = 0Aa> rec , V c0 = 0, r = 0, 
& eQ = 0.5, & es = 0.5, & gQ = 0, ^„ = and ® eQ = 0.7. 



The dipole-dipole potential V responsible for the non- 
linearity is switched on at time t = and the atoms are 
assumed to be left in their excited state after being 
diffracted by the grating. In that case, the conjugate 
beam builds up only in the electronic ground state as 
expected from the analytical solution in the undepleted 
regime [13]. The partial Schrodinger waves first 
undergo a transient, with discontinuities in their slope 
after the beams cross the condensate an integer number 
of times. 2 The kicks become less and less pronounced as 
a new steady state corresponding to V ^ is established. 

Figure 3 shows how the fraction of excited state atoms 
in the reflected beam influences the intensity of the 
conjugate beams. As follows both from analytical con- 
siderations [13] and from our numerical simulations, no 
conjugate beams are produced if the atoms exit the 
grating in their ground state. This is simply because in 
that case, there is no spontaneous emission, and no 
dipole-dipole interaction between the atoms. Note the 
asymmetry between the excited and ground state atomic 
conjugation, which results from the asymmetry in the 



pump beams, the condensate being assumed to consist 
of ground state atoms. Of course, the results for large 
fractions of excited atoms should be regarded with cau- 
tion, due to the neglect of the incoherent effects of 
spontaneous emission in these results. This is the object 
of the next section. 




0.25 
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This is an artifact of our numerical algorithm, which starts propagat- 
ing the states from the two boundaries at z = and z = L after V is 
switched on. 



Fig. 3. Strength of the conjugate beams|£$ e (-K, 0, t)\ 2 — dashed 
line, and | % g (-k, 0, t)\ 2 — solid line, in steady state, as a function of 
the fraction 17 of excited atoms in the reflected beam. k = Jc , k= 10k , 
L = 40/&0, Vo = 0.4w rec , V c0 = and r = 0. 



4. Spontaneous Emission 
4.1 Master equation approach 

In order to account for the effects of spontaneous 
emission, the Hartree nonlinear Schrodinger equation, 
Eq. (6) must be replaced by the nonlinear master 
equation, Eq. (1), as already mentioned. In analogy with 
the coherent case, we now expand the effective single 
particle density matrix element pw(r, r\ t) as 



p<? (r, r\ t) = p& (r, r\ t)c ik ^ + p# (r, r\ t)c~ ik ^ 
+ p f i b t (r, /•', f )e* (z+ ° + p'f (r, r\ f )e"* (z+0 . (10) 



This decomposition yields 16 coupled paraxial-like 
partial differential equations of the spatial variables 
(y,z), (y',z') and t. Their explicit form is given in 
Appendix A. These equations, being first-order in z and 
z\ require boundary conditions on the four sides of the 
square z, z' €. [0, L]. However, the boundary conditions 
Eq. (9) on the wave functions in the coherent case yield 
boundary conditions for the density matrix at the cor- 
ners of this square only. The additional boundary equa- 
tions required for the sides of the square can be derived 
from the Hartree nonlinear Schrodinger equation with 
the complex potential V + V c , as shown in Appendix A. 
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We have numerically solved the Hartree-Fock nonlin- 
ear master equation in two dimensions for various spon- 
taneous decay rates and fractions r\ of excited atoms at 
the exit of the beam splitter. As we are interested in 
particular transverse directions of the wave propagation 
we introduce a spatial Fourier transform of the atomic 
density operator in the transverse dimension 

p$(k y ,z;k f y ,z\t) = 

Idkydk'y p$ (y 9 z;y r , z\ t)e lk y y e*> y ', (1 1) 




0.25 



0.75 



where a and (3 stand for the forward and backward 
propagating wave. 3 

Figure 4 shows the probabilities p% (-k, 0; -k, 0, t) 
and pn (0, 0; 0, 0, t) as a function of time, for the case 
where the atoms exit the diffraction grating in their 
excited state. As such, it is the same as Fig. 2, except that 
now the spontaneous emission rate has been set equal to 
r= 0.5 and the imaginary potential is V c = 0. 1, in recoil 
units. Clearly, the amount of phase conjugation is sub- 
stantially reduced as compared to the coherent situation. 
This is further illustrated in Fig. 5, which shows the 
steady state phase conjugate probability as a function of 
the spontaneous decay rate. 

Already for a decay rate equal to the recoil frequency, 
phase conjugation has all but disappeared. 4 
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Fig. 4. Master equation simulations in the incoherent regime. Time- 
dependence of the backward propagating atomic beams 
p g ^(0, 0; 0, 0, t) — solid line, pf g (-k, 0; -k, 0; t)— dashed line, for 
K = k , k= lOifco, L = 40/k Q , V = 0.2a> rec , V cQ = 0.1 a> rec , r=0.5w rec , 
Vnt = *o, &eo = 0.5, & es = 0.5, Pgo = 0, & gs = and 5^=0.7. 



The parameters of our simulations were largely dictated by com- 
puter limitations. It is expected that phase conjugation will survive 
larger spontaneous emission rates if larger dipole-dipole potentials and 
longer propagation distances are considered. 

As it should be, our numerical results agree with the coherent ones 

for r=v c = o. 



Fig. 5. Influence of r on the strength of the conjugate beam 
p^(-ic, 0;-k, 0; t) in steady state K = k , k=\0k , L = 40/& , 

, nt = k Ql & eQ = 0.5, ^ = 0.5, ^ g0 = 0, 



V = 0.2w rec , V c0 
& gs = and 5^o = 



= 0. \a re 
0.7. 



4.2 Hartree Monte Carlo Simulations 

The numerical solution of the nonlinear master equa- 
tion leads to substantial computer memory require- 
ments. In a recent paper [24], we have proposed an 
alternative approach to its brute force solution, based on 
the average over "quantum trajectories" resulting from 
the solution of a Hartree nonlinear Schrodinger equa- 
tion, interrupted by quantum jumps. The validity of this 
approach, which is inspired by the now well-established 
Monte Carlo wave functions technique [25-27], can 
however only be proven provided that the Hartree ap- 
proximation is exact, which holds strictly only for con- 
densates. In addition, the two-body contribution to dis- 
sipation must be of the form of a complex potential, 
which is the case in the present situation. A direct com- 
parison of the Hartree Monte Carlo and master equation 
results should therefore give us not only an indicator of 
the practical validity of the former method, but also a 
sense of how spontaneous emission destroys the conden- 
sate. 

In order to implement this approach, we generalize 
the Hartree nonlinear Schrodinger equation, Eq. (6) to 
include the imaginary potential V c , so that it becomes 



a d£<g} = [ d2 ( a \ H \2)<l>(2) 



>/■ 



+ Jdld2d3 (a, l\V+ V c \2, 3) <f>*(l) <f> (2) <f> (3) . 

(12) 

Expressing furthermore the linear Liouvillian L nc in 
Eq. (1) as the Lindblad form 
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L nc p (t) = -- 2 ( Cl C m p(t) + p (t) Cl C n 



+ 2 c mP (t)d 



(13) 



where the C m 's are system operators appearing in its 
coupling to the continuum of electro-magnetic field 
modes, the Hartree Monte Carlo algorithm proceeds as 
follows [24]. We first formally reexpress the nonlinear 
Schrodinger equation, Eq. (12) as 



In case no jump occurs, the state vector evolves under 
the influence of the nonlinear and nonhermitian effec- 
tive Hamiltonian i/ eff , so that 



m (f+8f)>~ \<j>(t))- — ht 



h 



and the new wave function is 



(18) 



m\cj>(t)) = H nl [\4,(t))]\4,(t)). 



(14) 



| <j>{t + ht)) = | <j> m (t + 8f))/Vl-8/7 . 



(19) 



where the nonlinear "Hamiltonian" H ni [\<j> (t))] 
is a functional of the effective single-particle state 
vector \<j>(t)), see Eq. (12). We then decompose the 
Liouvillian L nc into a contribution H d added to H ni to 
yield the effective nonhermitian "Hamiltonian" 



H c v[\<l>)]=H nl [{</>)] +H d 



= H nl [\<j>)]- l j- 2 Cl C m 



(15) 



and a linear "fill-up term," see Eq. (13). As in the 
usual case, H^ is responsible for the Schrodinger-like 
evolution of the system, while the "fill-up term" causes 
randomly distributed quantum jumps. 

The probability 8/? for one of these jumps to occur is 
given for the small time interval ht by 5 



g _iht_ (4>(t)\H d -Hl\$(t)) 

P ~ h <0(O|^(O> 



= 8r 



<<«ol2 elegit)) 

m 

m)\cf>(t)) 



; 2 &p» 



(16) 



We performed Monte Carlo simulations within the 
slowly varying envelope approximation introduced in 
Sec. 2 and for the same beam configuration as in the 
preceeding example (see Appendix B). As pointed out in 
Ref. [24], the accuracy of the Monte Carlo simulations 
as compared to the Hartree-Fockmaster equation 
description is expected to depend on the amount of 
dissipation in the problem. Figure 6 shows the relative 
difference between the two methods as a function of the 
spontaneous emission rate r. These results confirm that 
the Hartree Monte Carlo simulations work surprisingly 
well, although the discrepancy between the two 
methods increases with 7", as expected. 




A jump is said to occur if a quasi-random number e, 
uniformly distributed between and 1, is less than 8/?. 
The post-jump wave function is chosen amongst the 
various possible final states C m \(f)(t)) according to the 
probability law 8/? m /8/?, with 



|(/>(r + 80) = C m |^(0)/(8/ ?m /80 1 



(17) 



Note that we keep track of the norm of the state vector explicitly, 
since only a partial volume of the system is considered explicitly, and 
thus the flux, rather than the norm, is conserved in our numerical 
solution. 



Fig. 6. Comparison of the Monte Carlo and master equation simula- 
tions. Influence of Fon A, the standard deviation between the steady 
state results of the Monte Carlo and master equation predictions, 
normalized by the Master equation result. 



We attribute this surprisingly good agreement to the 
fact that the type of dissipation that we are considering 
here is however quite special, in that the only effect of 
the nonlinear Liouvillian L ni is in the form of an imagi- 
nary two-body potential V c . It corresponds to a modifi- 
cation of the spontaneous decay rate resulting from the 
interference of the spontaneous emission probability 
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amplitudes due to the undistinguishability of the two 
particles. In other words, it describes a collective decay 
mechanism (albeit as a result of the Born-Markov 
approximation, only two-atom effects are retained). Due 
to this collective nature, it is reasonable to expect that 
this imaginary potential might not do too much violence 
to the Hartree ansatz. 

The only other incoherent contribution to the nonlin- 
ear master equationis a single-particle Liouvillian L r . 
But due to the indistinguishability of the particles in the 
sample, it describes a decay mechanism that occurs for 
all N atoms simultaneously. As a result, it leaves all 
atoms in the same state, albeit now described by a 
product of N identical single-particle density operators. 
It appears then that for the very specific form of decay 
we are considering, the Hartree ansatz might not be as 
bad as otherwise expected, a conjecture largely con- 
firmed by the result of Fig. 6. 



5. Summary and Conclusions 

In this paper we have discussed the possibility of 
observing atomic phase conjugation off Bose conden- 
sates, and to use it as a diagnostic tool to determine the 
spatial coherence properties and lifetime of that conden- 
sate. Since phase conjugation results from the scattering 
of a partial matter wave off the spatial grating produced 
by two other waves, it offers a very natural way to 
directly measure these coherence properties, and as 
such provides an attractive alternative to the optical 
methods proposed in the past. 

As expected, incoherent mechanisms such as sponta- 
neous emission rapidly destroy phase conjugation, 
thereby demonstrating that the proposed technique 
should easily distinguish a condensate from an incoher- 
ent sample of the same density. We have also compared 
Hartree-Fock master equation results to those based on 
Hartree Monte Carlo wave function simulations, thereby 
obtaining a measure of the destruction of the condensate 
as a function of the rate of spontaneous emission. 

We conclude by noting that in replacing the nonlocal 
dipole-dipole potential by a local interaction, the grat- 
ings established in the condensate are probed only 
locally. It will be interesting to investigate if more infor- 
mation can be gained by treating V as a nonlocal poten- 
tial. In particular, one can expect that effects such as 
spatial diffusion in the sample should be accessible by 
such nonlocal sampling. 



6. Appendix A. Master Equation 
in the Slowly Varying Envelope 
Approximation 

In analogy to the wave function description in terms 
of forward and backward propagating beams each of the 
four density matrix elements p ee , p gg , p eg , p ge can be 
decomposed into four parts 

per (r, r\ t) = p& (r, r\ t)e ik(z ~ z,) + p$ (r , r\ t)* 4 *™ 

+ p& (r, r\ 0e^ +o + p# (r, r\ Oe"*^ , (20) 

where 

p& = F f (r,t)B*(r' 9 t). (21) 

The 16 coupled paraxial-like partial differential equa- 
tions obtained from the nonlinear master equation, 
Eq. (1) in the slowly varying envelope approximation are 
of the form 

i(Vl-V?) P i£(r,r';t) 

- iV* [(pi{(r, r; t) + p** (r, r, t)) pg(r, r'; t) + 

p%{r,r;t)p b J{r,r';t)} 

- iV* l(p&(r, r, t) + p** (r, r; t)) P S(r, r'; t) + 

pg(r,r',t)p*(r,r';t)] 

+ iVJHpg(r; r'; t) + p** (r\ r't)) pg(r, r'; t) + 
pV(r',r';t)p%(r,r';t)] 

+ iV i0 [{pg (r\ r'; t) + p»» (r\ r'; t)) pg{r, r'; t) + 
p b ((r',r';t)pg(r,r';t)l 



+ {er'\L nc pf e {t)\er), 



(22) 



where V, = V + V c0 . Similar equations can be obtained 
for all other components. 

These first-order differential equations in z and z' 
require boundary conditions on the four sides of the 
square z,z'=[0,L]. However, the boundary conditions 
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for the wave functions in the coherent case Eq. (9) yield 
boundary conditions for the density matrix only at the 
corners of the square. Ignoring spontaneous emission at 
the boundaries, we can express the density matrix 
elements as pg (y 9 z, y\ 0, t) = F e (y 9 z, t)F*(y\ 0, t), 
etc. In this expression, F e (y,z,t) undergoes a 
Schrodinger evolution, while F*(y\ 0, t) is a boundary 
value. Thus, the dynamics of the system on the sides 
of the square is obtained by considering a boundary 
equation for the master equation, which involves only 
the evolution corresponding to the change of 4 'half ' of 
the density matrix element. 

This boundary equation can be derived from the 
Hartree nonlinear Schrodinger equation with the poten- 
tial V + V c . In particular, we find for p{{ (r, r' b0 und; 



resembles atomic recoil due to the spontaneous emission 
event [26]. The effective Hamiltonian // eff [\<fi )] is then 



H eff [\cf>)]=HA\4>)]- l -^fo- + o- 



and yields the quantum jump probability 



3jy l<MQl 2 

I*(0| 2 



8/7 = ^-8; 



£ 



d*d-* 2 ) 



2 §Px- 



(25) 



(26) 



( — + 2k — j p{{ (r, r' bound ; t) = 

i Vlp2(r,r bound; t) 

*V% HP* (T,r\t) + p e b g b (r, r ; 0) P$ (r , r' bound ; + 
p^ (r , r ; t) p h g { (r, r' boU n d ; 01 

W% lipfg ( r > r \ + p// ( r ' r ; 0) p// (r, /" bound ; i) + 



?// 



p/| (r, r ; p*/ (r, r' bound ; 01 



(23) 



Finally, the post-jump wave function turns out to be 



e-'>*V- 



|0(f + 8f)> = - 



10(O> 






«)i 2 



l*(0l 2 



(27) 



i.e., <fi e (t + 8?) = and 

<M* + 80 = e- ! ?^Vi^(0)/V|^(0| 2 /|</>(0| 2 . 



Here r' bound is the position coordinate on the boundary, 
** bound = (x\y',0) for the forward propagating compo- 
nent and r' bound = (jc', y',L) for the backward 
propagating component. Note that this situation is 
different from the case of confined systems, where 
pet (y 9 z\y\ z'bound, t) = on the boundary. 



7. Appendix B. Explicit Expressions 
for the Hartree Monte Carlo 
Simulations 



The explicit form of the incoherent Liouvillian L n 
Eq. (1) leads to the system operators 



in 




(24) 



where 3(1- ^f 2 )/4 is the dipole radiation pattern for the 
two-level atoms and ;^=cos#, being the angle be- 
tween the quantization axis and the direction in which the 



spontaneous photon is emitted and exp 



±ife(k)r 



In case of decay the wave functions are given by 
Eq. (19). 

These simulations are carried out within the slowly 
varying envelope approximation introduced in Sec. 2. 
The expressions for the slowly varying amplitudes 
B e (r 9 t+bt\B g (r 9 t+bt\F e (r,t+bt) and F g (r,t+bt) can 
be found easily from the Eqs. (17) and (19). 
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